流体力学与物理导引神经网络
摘要
现代科学发展至今，众多复杂的关键应用依赖着对流体运动的精确预测：大到全球气候的预测，小到飞机机翼与发动机的设计，上至航天火箭液体燃料对本体造成的稳定性问题，下至海浪能源的提取等。然而，由于流体力学方程数值求解的复杂性，基于传统数值方法的高分辨率长周期数值模拟运算量大且难以保证数值稳定性，流体力学问题的数值模拟因此成为应用科学领域中的重大挑战之一。随着数字时代带来的冲击，新方法、新尝试不断更新迭代，利用深度学习解决复杂的数学问题已成大趋势。在上一期中，我们已经大致介绍了近些年研究热点物理导引神经网络(PINN：Physics Informed Neural Network)的基本研究与应用。本期将更进一步，深入介绍PINN方法在流体计算中的应用，以及该方法与传统数值方法的差别。
感谢彭峻博士对本期流体力学部分的悉心建议与修改。
得益于硬件与软件框架的高速发展，神经网络在不同领域上获得了许多引人注目的突破：从最初的全连接网络，到专用于自然语言处理和时间序列的RNN[1]，LSTM[2]，Transformer[3]等网络结构；有运用于计算机视觉领域里基于卷积网络(Convolutional Neural Network: CNN)的Inception[4]，ResNet[5]，EfficientNet[6]等衍生结构；还有专门精于从氨基酸系列就能精准预测蛋白质3D结构的AlphaFold[7]。但在大部分应用中，这些神经网络多是基于现有的数据集进行监督学习。例如使用标记好狗与猫的图片来进行动物种类分类的训练。然而数据的收集与标记往往需要大量的时间与资源，在一些应用场合很难有效的达成。但是如果能够利用上一些先验知识，如依据物理规律或者其他原理，无需观测数据即可训练神经网络，那么将极大地提升深度学习的可应用范围。
在基于先验知识或原理的神经网络中，最有潜力的方法当属Physics-Informed Neural Network (PINN)[8]。该方法在损失函数里加入了基于物理方程的约束条件，使其能够有效地利用物理先验知识，甚至无需使用任何观测数据也能进行相关的物理现象的预测。正如上一期提到的，PINN方法已经在物理、化学，生物等众多领域展示了应用的潜力。其中一个受人瞩目的应用领域当属流体力学，例如其已被应用于包括流体图像的超分辨率问题[9]。
由于流体运动的高度非线性，人们对流体的诸多复杂现象，如湍流，泡沫形成与破灭等知之甚少，所以有关流体长期运动的计算一直都没有好的解决方案。所以就算在深度学习盛行的如今，关于流体的传统数值计算方法的研究仍然相当活跃。本期将在简要介绍流体力学数学描述的基础上，对传统计算流体力学方法和基于神经网络的方法进行比较。
流体一般可由两种方式进行描述：以固定网格坐标来测量流体变量的欧拉 (Eulerian) 视角（图1蓝色）和以追踪流体粒子运动轨迹而随时间变化的网格为拉格朗日(Lagrangian) 视角（图1红色）。在二维空间里，在拉格朗日视角下的自变量为粒子本身的标签或参考坐标 a, b和时间τ。所求的因变量为粒子位置的横坐标x(a,b,τ)与纵坐标y(a,b,τ)，流体密度ρ(a,b,τ) 和压力P(a,b,τ)。这是最自然的流体表达方式，因为我们能追踪每一个粒子最初状态以及之后的演化，继而在实际数值计算中会拥有更好的能量守恒等特性。但随粒子的物理坐标跟着时间变化，可直接测量的坐标位置会受限制，可能需使用插值法来延伸粒子的物理位置，然而插值法会大大降低计算的精确度。因此拉格朗日视角更适合应用到无需计算流体变量在某些固定位置的问题中，例如液体与容器的耦合受力和高楼防震防风的阻尼器等问题上。
在使用欧拉视角来描述流体运动时，自变量为粒子位置的横坐标x、纵坐标y和时间t。所求的因变量为流体横向速度u(x,y,t)、纵向速度v(x,y,t)，流体密度ρ(x,y,t)和压力P(x,y,t)。欧拉描述适合于需要观测给定位置流体信息的问题，如对某地的风速与大气压的观察。著名的纳维尔-斯托克斯方程(Navier-Stokes equations)便是基于欧拉视角来描述流体的运动。
根据所考虑的问题不同，还可将流体分为可压缩与不可压缩两大类。在大部分情况下，流体的密度在空间里为均匀并且不随时间变化 (例如水)，即不可压缩流(incompressible flow)。在这种情况下，ρ(x,y,t) = ρ为已知常数，二维的纳维尔-斯托克斯方程可表示为：
这三个方程分别对应了流体横向动量，纵向动量和质量守恒方程。该系统的输入为x, y, t，输出为u, v, P。方程里的系数μ为动力粘性系数(dynamic viscosity)，g为重力加速度，ρ为密度以及已知形式的外力F(x,y,t)和G(x,y,t)。
由于动量方程中非线性项（也称之为对流项）的存在，除了在特定条件与假设下，该方程难以解析求解。需要强调的一点是，在不可压缩流里，流体压力P(x,y,t) 本身的绝对数值是毫无意义的，因为方程中只有压力梯度的存在，这会造成在数值计算中的一些困难。首先，流体压力是没有时间导数的项，所以这会造成P(x,y,t) 的时间离散化无法被写成传统数值算法要求的形式来进行时间积分。也因此，在解流体质量守恒的方程中通常需要解超大线性方程组，而且这些线性方程组在大部分情况下是奇异 (singular) 矩阵，会带来巨大的计算成本。另外，空间上的中心差分离散会造成流体压力的奇偶失联 (图2)从而产生一种叫‘棋盘’(checkerboard)的非物理现象，它会造成迭代过程中误差的不断扩大，因此一般会需要使用交错网格(staggered grid)或拟压缩法(artificial compressibility method)等方法来解决这类问题。
由克雷数学研究所 (Clay Mathematics Institute)提出的每道价值100万美金的世界七大数学难题之一就有证明三维的不可压缩纳维尔-斯托克斯方程的解的存在与平滑性[10]。至今这七大难题只有一题，庞加莱猜想 (Poincaré conjecture)被俄罗斯数学家格里戈里·佩雷尔曼 (Grigori Perelman) 所破解。他本人并没有因攻破了这道著名的数学难题后接受奖金，而是出于对数学的喜爱一心想证明自身感兴趣的数学难题。
目前，二维中该方程组解的存在与平滑性已经被数学家证明。这是利用了二维空间里流函数 (stream function) 的标量性，让u=∂ψ/∂y和v=-∂ψ/∂x，ψ(x,y,t)则是标量流函数。使用流函数的好处在于它直接满足了质量守恒的条件并且把二维的速度场降为一维。这使得纳维尔-斯托克斯方程转换成了仅仅有一项的标量方程
像纳维尔-斯托克斯方程等偏微分方程(PDE)是以连续变量坐标来定义的。例如，即使横坐标x的取值为[0,1]之间的实数，也会存在无限多的数值量，无法在计算机上完整表述。因此需要将其离散化成有限的数值点，例如将x分割成{0, 0.1, 0.2, 0.3, … , 1.0}的离散网格点。
传统的流体力学数值模拟方法，可大体上分为两种，基于网格(mesh-based)和无网格(mesh-free)。其中网格式方法适用于欧拉视角的计算，而无网格方法因流体粒子位置可能会使网格严重扭曲，所以更多用于拉格朗日视角的计算。一般情况下，网格式方法可以更准确地进分辨流动的多尺度现象，例如湍流；而无网格方法可以更好地计算出粒子的位置信息，如气液界面、激波等。
常用的网格算法有基于在网格上做差分来求数值导的有限差分法(finite difference method)。还有更进阶且基于数学理论的有限元法(finite element method)，将整个空间域分解成用简单的基函数来代表的计算域。如果基函数被设为定值全域的傅里叶基时也称之为谱方法(spectral method)。当基函数为分段常数函数时，该方法则是一直以来在计算流体动力学领域里最主流的数值算法，有限体积法(finite volume method)，因其数值方案从原理上保证了流体质量的守恒。当流体的表面积远远小于体积时，通常会用基于表面积的算法，例如边界元法(boundary element method)，将大部分计算转移到表面积上而减少总体的计算量。
无网格方法中最经典的莫过于光滑粒子法(Smoothed-Particle Hydrodynamics: SPH)，在拉格朗日视角上描述流体，并基于流体粒子之间的相互影响来计算每个粒子的位置、速度以及其他变量。
此外，还有一些结合了网格法和无网格方法特点的方法。例如隐式粒子法(Fluid-Implicit-Particle: FLIP)，在拉格朗日视角上描述流体粒子运动，而粒子之间互动的计算是在网格上实现的。又如任意拉格朗日-欧拉方法(Arbitrary Lagrangian-Eulerian: ALE)是一种基于有限体积/有限元网格的方法，适合于对界面追踪精度要求较高的问题。
最近，在一篇备受关注的文章中[11]，作者针对模拟大范围海洋运动问题提出了一种融合了用于计算流体粒子细节运动的隐式粒子法和用于计算大范围海面的边界元法的新方法，用于对快艇在海洋中穿梭的高效准确模拟。
应用在流体计算的深度学习方法其实有很多：包括上上期提到使用图神经网络(GNN)观测数据来进行流体的模拟，还有上期文章里介绍的物理导引神经网络(PINN)。相较而言，因PINN方法对描述物理的方程无需使用观测数据而能直接拟合的优势，使其成为热门的计算流体的无网格方法。相对于传统数值方法的迭代算法，PINN方法使用神经网络直接对整个空间与时间域进行拟合。这里先大致回顾一下PINN方法。
PINN方法的核心思想很简单，主要就是将描述物理或某些其他过程的方程加进损失函数(loss function)里，因此可以写为
在这里，x是神经网络的输入而u则是输出。例如在三维流体中，x=(x,y,z)和 u=(u,v,w)对应了三维的空间坐标与速度；R则包含了需要解的方程或方程组。在流体力学里，R对应的可以是三维的纳维尔-斯托克斯方程，因此能写成
R1–R3是每个维度对应的流体动量方程，Rm则是质量守恒的方程。不过我们还需要初始与边界条件(例如：流体不可穿墙)来使整个问题适定(well-posed)从而保证系统存在唯一解(unique solution)。因此整个损失函数表达为
其中Q和S各自包含了初始条件与所有边界条件。例如当U为初始速度且b为边界值时，Q=∑u(x,t=0)-U(x)和S=∑u(x=0,t)-b。另外α和β这两个权重系数的调整能帮助优化损失函数的收敛。该损失函数的优化则一般会使用一阶的ADAM优化器或者拟二阶的L-BFGS-B优化器来获得更快的收敛效果。实际应用里，α和β的值会选为大于1来加强对初始与边界的拟合从而更好地拟合内域的方程，不过也存在根据Q和S的值与梯度而进行动态更新的方法。
总体来说，Q和S项对应的是数据拟合(data fit)而R项对应的是神经网络的正则化(regularisation)。也因此会造成损失函数的梯度下降不均匀和变化过快(刚性) ，进而影响整个优化过程的收敛并导致最终效果不佳。为了解决这个问题，一种方法是在方程模型里强制加入满足初始与边界条件的函数。例如设一个新的量
在这里，u是神经网络估算的输出，B是边界或初始条件的值，而D则是一项在初始时间或者边界位置时值为0的函数，因此η能完美地满足初始与边界条件。因此如果将η 代入到原有的方程里，则会产成一组新的方程让神经网络拟合。拟合出的结果经过以上的公式转换后则是最终的解。其实这里核心的思想是与其让优化器去寻找最优的解来满足初始与边界条件，不如强制把这两个条件加入在方程里来缩小优化器的搜索空间。
另外，由于流体计算的难度大并且成本很高，PINN方法对长时间的流体运动模拟通常是无法完成的，因为损失函数的优化无法得到收敛。也因此产生了迭代型的PINN方法，在时间维度上进行切割分块，然后把前一时间段算出的解加入到当前时间段的损失函数里进行迭代拟合而使PINN方法向后兼容(bc-PINN)[12]。
最后我们来讲讲PINN方法使用的网络结构。最传统的结构有全连接网络(图4)，该网络结构相对简单，并且可以很方便地实现，所以这也是目前PINN方法使用最多的网络结构。另外还有利用了在transformer上类似于attention机制的网络(图5)，因该网络结构能更好地保证损失函数在训练过程中梯度的一致，来获得更好的收敛速度与效果。除此之外还有基于CNN网络结构的，将方程从物理空间转换到规则的参考空间里(类似于图片像素)，从而使CNN能够应用其中(图6)。函数导数由在参考坐标上的有限差分法来实现，因此边界条件可以使用填充来直接满足。不过目前最适合PINN方法的网络结构还没有统一的意见，并且不同网络结构带来的效果都与相关问题有着密切的关联。
相较于传统的数值方法，基于神经网络的PINN方法在流体应用中的优点在于，计算上无需一步步迭代便可在整个时间领域内对方程进行拟合，并可以在短时间内更精准高效地捕捉到流体的非线性现象。使用上也更加简易，当需要添加或者更改方程里的某项或边界条件时，传统方法因用了复杂的离散处理而不易更改，并且有可能会破坏已有数值方案的一些良好属性，但是PINN方法只需要更改方程的形式即可。
不过PINN因为是相对较新的方法，所以有许多地方还不够成熟，特别是在数值与收敛分析上还很欠缺。此外，虽然前文介绍了几个常见PINN方法使用的神经网络结构，但专门针对PINN方法的网络结构还没有被人们挖掘。
针对这些不足点，我们认为若能利用到系统对称性来改造网络结构，可能会让PINN方法在未来获得一些突破。对称性在数学上是一个相当重要的概念，它能帮助把非常复杂困难的问题化简。在二战中，英国著名的数学家图灵(Alan Turing)就是利用了当时德军使用的加密法编码里可逆的对称性，从而破解了德军的恩尼格码加密系统(Enigma machine)，使盟军获得在情报上的优势。卷积网络在计算机视觉的成功也是基于它利用了图像里物体的平移(translational)对称性。许多物理系统都存在对称性，例如在流体力学中，粒子拥有可重新标记的对称性，因为流体粒子的所有变量都不会因为改变标记系统而改变。除此之外，能在网络架构上保证流体质量守恒，例如类似传统有限体积法的神经网络，也可能会在计算流体运动上有重大帮助。
[1] D. E. Rumelhart, G. E. Hinton, and R. J. Williams, “Learning representations by back-propagating errors,” Nature, vol. 323, no. 6088, pp. 533–536, Oct. 1986, doi: 10.1038/323533a0.
[2] S. Hochreiter and J. Schmidhuber, “Long Short-Term Memory,” Neural Computation, vol. 9, no. 8, pp. 1735–1780, Nov. 1997, doi: 10.1162/neco.1997.9.8.1735.
[3] A. Vaswani et al., “Attention Is All You Need,” Jun. 2017, [Online]. Available: http://arxiv.org/abs/1706.03762
[4] C. Szegedy et al., “Going deeper with convolutions,” Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition, vol. 07-12-June, pp. 1–9, 2015, doi: 10.1109/CVPR.2015.7298594.
[5] K. He, X. Zhang, S. Ren, and J. Sun, “Deep Residual Learning for Image Recognition,” Dec. 2015, [Online]. Available: http://arxiv.org/abs/1512.03385
[6] M. Tan and Q. v. Le, “EfficientNet: Rethinking Model Scaling for Convolutional Neural Networks,” May 2019, [Online]. Available: http://arxiv.org/abs/1905.11946
[7] J. Jumper et al., “Highly accurate protein structure prediction with AlphaFold,” Nature, vol. 596, no. 7873, pp. 583–589, Aug. 2021, doi: 10.1038/s41586-021-03819-2.
[8] M. Raissi, P. Perdikaris, and G. E. Karniadakis, “Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations,” Journal of Computational Physics, vol. 378, pp. 686–707, Feb. 2019, doi: 10.1016/j.jcp.2018.10.045.
[9] H. Gao, L. Sun, and J.-X. Wang, “Super-resolution and denoising of fluid flow using physics-informed convolutional neural networks without high-resolution labels,” Nov. 2020, doi: 10.1063/5.0054312.
[10] C. L. Fefferman, “Existence and smoothness of the Navier-Stokes equation,” The millennium prize problems, vol. 57, no. 67, p. 22, 2006.
[11] L. Huang, Z. Qu, X. Tan, X. Zhang, D. L. Michels, and C. Jiang, “Ships, Splashes, and Waves on a Vast Ocean,” Aug. 2021.
[12] R. Mattey and S. Ghosh, “A Physics Informed Neural Network for Time-Dependent Nonlinear and Higher Order Partial Differential Equations,” Jun. 2021, [Online]. Available: http://arxiv.org/abs/2106.07606
[13] X. Jin, S. Cai, H. Li, and G. E. Karniadakis, “NSFnets (Navier-Stokes flow nets): Physics-informed neural networks for the incompressible Navier-Stokes equations,” Journal of Computational Physics, vol. 426, Feb. 2021, doi: 10.1016/j.jcp.2020.109951.
[14 S. Wang, Y. Teng, and P. Perdikaris, “Understanding and mitigating gradient pathologies in physics-informed neural networks,” Jan. 2020.
[15] H. Gao, L. Sun, and J.-X. Wang, “PhyGeoNet: Physics-informed geometry-adaptive convolutional neural networks for solving parameterized steady-state PDEs on irregular domain,” Journal of Computational Physics, vol. 428, Mar. 2021, doi: 10.1016/j.jcp.2020.110079.
往期推荐
壁仞科技研究院作为壁仞科技的前沿研究部门，旨在研究新型智能计算系统的关键技术，重点关注新型架构，先进编译技术和设计方法学，并将逐渐拓展研究方向，探索未来智能系统的各种可能。壁仞科技研究院秉持开放的原则，将积极投入各类产学研合作并参与开源社区的建设，为相关领域的技术进步做出自己的贡献。